Dynamic phase separation of fluid membranes with rigid inclusions 
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Membrane shape fluctuations induce attractive interactions between rigid inclusions. Previous 
analytical studies showed that the fluctuation-induced pair interactions are rather small compared 
to thermal energies, but also that multi-body interactions cannot be neglected. In this article, it is 
shown numerically that shape fluctuations indeed lead to the dynamic separation of the membrane 
into phases with different inclusion concentrations. The tendency of lateral phase separation strongly 
increases with the inclusion size. Large inclusions aggregate at very small inclusion concentrations 
and for relatively small values of the inclusions' elastic modulus. 

PACS numbers: 87.16.Dg, 64.75. 4-g, 05.10.Ln 



I. INTRODUCTION 

Biological and biomimetic membranes consist of a lipid 
bilayer with various types of macromolecules such as pro- 
teins 1^, ||. Many of these macromolecules are incor- 
porated in the bilayer, others are covalently bound or 
adsorbed to the membrane. The membranes are fluid 
and often tend to phase-separate and to form domains or 
'rafts' with different molecular composition. In biologi- 
cal membranes, the presence of domains has been linked 
to specific functions in signaling budding [Q, or cell 
adhesion 1^, ^. In some cases, the domain formation is 
caused by a separation of the lipid bilayer into phases 
with different lipid composition |0, ||. In other cases, the 
phase separation of the membrane appears to be driven 
by attractive interactions between membrane inclusions 



Besides direct interactions such as van der Waals or 
electrostatic forces, membrane inclusions are subject to 
indirect interactions which are mediated by the mem- 
brane. Some of these membrane-mediated interactions 
are static, i.e. they arise from local perturbations of the 
bilayer structure or shape around the inclusions. Trans- 
membrane proteins which exhibit a hydrophobic mis- 
match with respect to the lipid bilayer cause a pertur- 
bation of the bilayer thickness. This thickness pertur- 
bation has been found to induce attractive interactions 
between the proteins |ll|, ^ |l^ . Similar interactions 
due to membrane thickness perturbations have also been 
proposed for adsorbed particles Membrane inclu- 

sions with conical shape 16j_17,_^, 19, ^ or 
membrane-anchored polymers ]227p3| cause local per- 
turbations of the membrane curvature which induces at- 
tractive or repulsive interactions. 

Other indirect interactions are dynamic, i.e. they are 
mediated by shape fluctuations of the membrane. In this 
article, rigid membrane inclusions are considered which 
interact due to the suppression of membrane shape flue- 
tuations H |4[ H, 1^ |8[ H, 0, see Fig. |. 



Fluctuation-induced interactions have also been found 
for specific receptors or stickers which locally bind op- 
posing membranes during adhesion ||3l| , |3^ , |3^ . 

The fluctuation-induced pair interactions of rigid 
membrane inclusions have been studied intensively |15| , 
|0|, H, H, H, 1^, |0). The pair interaction of rigid 
disks with radius r and distance L was found to be 



G{L) = -6kBT{r/L)^ + 0[{r/Lf] 



(1) 



in the absence of an external membrane potential and 
lateral tension |]l5[ Here, 0[{r/L)^] stands for terms 
of sixth or higher order in r/L. The prefactor of the 
leading term in Eg. p) was first 12, and later corrected to 
6 , see also |l^, |of |3^ . Since the distance L has to 
be larger than the inclusion diameter 2r, this term is only 
a fraction of the thermal energy ksT. However, higher- 
order terms in (r/L) might be relevant at small inclusion 
separations and thus contribute to the phase behavior of 
a membrane with many inclusions. Deducing the phase 
behavior of such a membrane is also complicated by non- 
trivial multi-body interactions pTj . 

In this article, it is shown numerically that membrane- 
shape fluctuations indeed lead to the aggregation of rigid 
inclusions. The phase behavior of a discretized mem- 
brane with rigid inclusions is considered in Monte Carlo 
simulations. The membrane consists of quadratic patches 
with linear extension a which corresponds to the smallest 
possible wavelength for bending deformations. Computer 
simulations of molecular membrane models indicate that 
this length scale is about 6 nm for a lipid bilayer with 
a thickness of about 4 nm [Q. Above a critical value 
K* for the stiffness modulus of the inclusions, the mem- 
brane is found to separate into an inclusion-rich and an 
inclusion-poor phase. The aggregation tendency of the 
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FIG. 1: Rigid inclusion in the lipid bilayer 



FIG. 2: A membrane segment containing inclusions with the 
size of one lattice site. The segment consists of 5 x 4 mem- 
brane patches which are labeled by the lattice sites i. Mem- 
brane patches with black inclusions correspond to occupation 
numbers = 1, while grey patches without inclusion have 
occupation numbers rii = 0. The local deviation of the mem- 
brane from the white reference plane is denoted by li. 




inclusions strongly increases with the size Q, which is 
reflected in a decrease of the critical stiffness K* with 
the inclusion size. Large inclusions also aggregate al- 
ready at relatively small inclusion concentrations. Here, 
quadratic inclusions with a size Q of 2 x 2, 3 x 3, or 4 x 4 
membrane patches are considered, which extends previ- 
ous results on smaller inclusions with the size of a single 
patch For these inclusions, the critical stiffness is 

found to decrease according to the power law K* ^ Q'^ 
with exponent c ~ —0.70. 



II. GENERAL MODEL 

The configurations of a fluctuating membrane with in- 
clusions can be described by a field I for the membrane 
shape and a concentration field n for the inclusions . 
For a membrane which is on average planar, the mem- 
brane shape is usually given by the deviation l{x,y) from 
a reference plane with coordinates x and y. Here, we 
discretize the reference plane into a square lattice with 
lattice constant a which corresponds to the smallest pos- 
sible wavelength for bending deformations. The inclusion 
positions are then given by occupation numbers ti^ = 
or 1 where rii — I denotes the presence of an inclusion at 
the lattice site i of the reference plane, see Fig. ||. 

In the absence of inclusions, the discretized bending 
energy per lattice site can be written as 



(2) 



where Ko is the bending rigidity of the lipid bilayer, and 
hi£i,x + Ci^y) is the local mean curvature of the membrane 



561. Here, 



— ij'x+a.y ~i~ ^x — a,y '^^x,y)/^ 
{^x.y+a ^x.y — a '^^x,y^/^ 



(3) 
(4) 



are the discretized curvatures in x- and y-direction at the 
lattice site i with coordinates {x,y). The rigid inclusions 



FIG. 3: Top view of membrane segments with inclusions of 
size Q = 1, 2 X 2, 3 X 3, and 4x4, respectively. The inclusions 
are shown in black. 



here are characterized by the elastic energy per site 



(5) 



with the stiffness modulus K g^. For K oo, such 
inclusions are completely rigid and suppress any local 
curvature at the inclusi on p osition similar to the rigid 
disks or rods studied in |15|, ^ [30| . In contrast, in- 
clusions with increased bending rigidity as considered in 
Ip^ , |29|| only suppress fluctuations of the total cur- 
but not saddle-type fluctuations with 



vature Ci^x 

The grand-canonical Hamiltonian of a membrane con- 
taining inclusions with the size of one lattice site can be 
written as 

i 

+ V{k)} (6) 

were fj, is the relative chemical potential of the inclusions, 
and V{li) is the external membrane potential. On lattice 
sites with occupation numbers rii — 1 indicating the pres- 
ence of inclusions, the elastic energy is given by 5/ (/) . On 
lattice sites with = 0, the elastic energy is the energy 
of the lipid bilayer. 
The larger inclusions considered here have an area of 
Q = 2x2, 3x3, or4x4 membrane patches or lattice 
sites, see Fig. |[ The elastic energy at every lattice site 
of an inclusion is given by Eq. (|^). The larger inclusions 
thus can be seen as quadratic arrays of small inclusions 
with the size of one lattice site. The grand-canonical 
Hamiltonian for a membrane with larger inclusions can 
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be formally written as [B2[ 



-HQil, n}^Y.\ ^''(^) + E (^4(0 - ^.f (O) - M 
+ V{k)\+Y,W^,n,n, (7) 



where {il, . . . , iQ} denotes quadratic arrays of Q = 2 x 2, 
3 X 3, or 4 X 4 lattice sites. The position of an inclusion 
given by rii = 1 corresponds to one of the lattice sites 
occupied by the inclusion, e.g. the center of an inclusion 
with size Q = 3 x 3. The hard-square interaction 



Wij = oo for j in Af 
— otherwise 



(8) 



prevents any overlap of inclusions. Here, Af denotes the 
exclusion area of an individual inclusion with size Q at 
lattice site i. 

In the following, the external potential of the mem- 
brane is taken to be the harmonic potential 



V{k) = 



(9) 



with potential strength m. The harmonic potential intro- 
duces an additional length scale, the correlation length 
^ = (4a^Ko/m)^/'* for the deviation field /, see, e.g., 
Ref. [g^. Membrane fluctuations on length scales larger 
than the correlation length ^ are suppressed by the har- 
monic potential, while fluctuations on smaller scales are 
governed predominantly by the elastic energy of the 
membrane. 

The membrane model defined by the Eqs. (|^) to (|^) 
has four characteristic dimensionless parameters, as can 
be shown be introducing the rescaled deviation field 



z = {l/a)y/Ko/{kBT) 



(10) 



These parameters are the ratio K/ko of the inclusion 
modulus and the bare membrane rigidity, the dimension- 
less chemical potential fj,/{kBT) for the inclusions, the 
rescaled potential strength rh = ma'^/Ko, and the inclu- 
sion size Q. 



III. MONTE CARLO SIMULATIONS 

To deduce the phase behavior from Monte Carlo sim- 
ulations, the inclusion concentration X = Q{ni)/a'^ is 
determined as a function of the chemical potential /i for 
various values of the inclusion stiffness K, size Q, and the 
rescaled potential strength rh. A first-order phase tran- 
sition is reflected in a discontinuity of X{ii) at a certain 
chemical potential fitr- The two limiting values of X{p) 
at are the inclusion concentrations of the coexisting 
phases, an inclusion-rich and an inclusion-poor phase. To 
determine the inclusion concentration X at a given value 



of fj,, Monte Carlo simulations are performed with up to 
10^ Monte Carlo steps per lattice site on a lattice with 
120 X 120 sites and periodic boundary conditions. Each 
Monte Carlo step consists in attempted local moves of the 
rescaled deviation field z and of the concentration field 
n on all lattice sites. For rescaled potential strengths 
rh > 0.01 as considered here, the correlation length of 
the membrane is much smaller than the lateral extension 
of the lattice, and finite size effects are negligible. 

In Fig. ^ phase diagrams as a function of the inclusion 
modulus K are shown for inclusions with size Q = 1, 2x2, 
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FIG. 4: Phase diagrams for inclusions with size Q = 1, 
Q — 2 X 2, and Q = 3 x 3 as function of the inclusion con- 
centration X and stiffness K in units of the bare membrane 
rigidity Ko. The rescaled potential strength is m = 0.01. In- 
side the shaded coexistence regions, the membrane separates 
into an inclusion-rich and an inclusion-poor phase with con- 
centrations given by the lines of Monte Carlo points. The 
critical points are represented by stars. 
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FIG. 5: Critical inclusion stiffness K* as a function of the 
rescaled potential strength ifi for inclusions with size Q = 
2x2. The critical stiffness K* is given in units of the bare 
membrane rigidity Ko and increases with m since the harmonic 
potential suppresses membrane fluctuations. 



FIG. 6: Scaling plot for the critical inclusion stiffness K* as 
a function of the size Q in the limit of the rescaled potential 
strength m — > 0, see also Table 1. The straight line has the 
slope —0.70 which corresponds to the exponent c in Eq. (pT[). 



TABLE I: Critical values K* / k,o of the inclusions stiffness K 
in units of the bare membrane rigidity Ho for various inclusions 
sizes Q and rescaled potential strengths m. Values for K* / Ko 
in the limit m — > are obtained by extrapolation from the 
critical values at finite m shown in the columns 2 to 4. 





rh — 1 


m = 0.1 


m = 0.01 


m ^ 


Q 


= 1 


300 ± 10 


225 ± 10 


205 ± 10 


198 ± 10 


Q = 


2x2 


100 ±2 


81 ±2 


77 ±2 


76 ±2 


Q = 


3x3 


55 ±2 


44 ±2 


43 ±2 


42.9 ± 2 


Q = 


4x4 


37 ±2 


29 ±2 


28 ±2 


27.9 ± 2 



length C leads to finite size effects in Monte Carlo simu- 
lations with m = 0, the critical stiffness K* for m — > is 
determined here by extrapolation, see Table 1. For the 
inclusion sizes considered in this article, the functional 
dependence of the critical values K* for to ^ on the 
size Q can be approximated by the power law 

if*~Q= (11) 
with exponent c = —0.70 ± 0.01, see Fig. ^. 



and 3x3. The rescaled potential strength is to = 0.01. At 
points {X, K) inside the shaded 2-phase coexistence re- 
gions, the membrane separates into an inclusion- rich and 
an inclusion-poor phase. The inclusion concentrations 
in the coexisting phases are given by the lines of Monte 
Carlo data points, the critical points are represented by 
stars. The extent of the 2-phase regions strongly in- 
creases with the inclusion size Q, wich is reflected (i) 
in a strong decrease of the critical stiffness K* and (ii) 
in an increase of the width of the coexistence region with 
the inclusion size. For inclusions with size Q = 3x3, 
the membrane phase-separates already at very small in- 
clusion concentrations. 

The external harmonic potential suppresses mem- 
brane shape fluctuations on length scales larger than the 
correlation length ( = o(4/to)^/^. Since the phase sepa- 
ration of the membrane is driven by the fluctuations, the 
critical stiffness K* increases with the rescaled potential 
strength to, see Fig. |[ For weak potentials with small 
values of to, the critical stiffness is rather independent 
of TO and tends towards a limiting value, since the cor- 
relation length then is much larger than the average 
distance between neighboring inclusions. 

In the absence of an external potential, i.e. for to = 0, 
the critical stiffness K* of the inclusions only depends 
on the inclusion size Q. Since the increasing correlation 



IV. SCALING ANALYSIS 

Large inclusions with size Q > 1 can be seen as 
quadratic arrays of small inclusions with the size of one 
lattice site. These inclusions aggregate at much smaller 
values of the stiffness K than inclusions with size Q — I, 
see above. In order to understand this behavior, it is 
instructive to consider the free energy difference AF 
between the two membrane states: The uniform state 
in which the inclusions are more or less homogeneously 
distributed throughout the membrane, and the phase- 
separated state in which the inclusions are aggregated. 
If AF is negative, the membrane will be in its homoge- 
nous state; if AF is positive, the membrane will phase- 
separate. For a given area concentration X of the inclu- 
sions, the free energy difference between the uniform and 
the aggregated state can be written as 

AF{K, Q) = AF,„t(if, Q) - TA^„„(Q) (12) 

in the absence of an external potential, i.e. for m = 
0. Here, AFint{K,Q) is the difference in the dynamic 
interaction free energy of the inclusions, which is in- 
duced by the shape fluctuations of the membrane, and 
ASmix{K, Q) is the difference in the entropy of mixing. 

The second term in Eq. ( [l^ ) is dominated by the en- 
tropy of mixing in the homogeneous state. For small 
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area concentrations X = Qljiijja} of the inclusions, this 
entropy of mixing is proportional to the number of inclu- 
sions, and ^Sraix scales as 



A5,„„(Q) - 1/Q 



(13) 



Eq. (|T^) simply results from the fact that the number of 
inclusions is proportional to l/Q for given area concen- 
tration X of the inclusions. 

The term AFint{K,Q) is dominated by the interac- 
tion free energy of the inclusions in the aggregated state. 
For Q > 1, the inclusions are rather densely packed in 
this state with an area fraction a^X larger than 0.9 (see 
Fig. ^) , and have contact with neighboring inclusions al- 
most along the whole circumference of length If 
one assumes that the interaction free energy in the ag- 
gregated state is proportional to the ratio of the inclusion 
circumference 4-y/Q and the area Q, the scaling form of 
the interaction free energy difference can be estimated as 

AF,„,(if,Q)~/(if)/v/Q~i^'"/v/Q (14) 

presupposing powerlaw-form for f{K). According to 
this estimate, the interaction free energy AFint decreases 
with the inclusion size Q proportional to Q^^^^, follow- 
ing the decrease of the 'surface-to-area' ratio of the in- 
clusions. The critical stiffness K* obtained from AF — 
then scales as 



-l/(2m) 



(15) 



Comparing with the exponent c ~ —0.70 from the Monte 
Carlo simulations (see Fig. H), one obtains 



m = -l/(2c) ~ 0.71 



(16) 



Thus, for the stiffness values K and inclusion sizes Q con- 
sidered here, the dynamic interaction free energy AFint 
of the inclusions appears to increase proportional to 

^0.71_ 



V. DISCUSSION 

In summary, I have considered the dynamic phase be- 
havior of a discretized membrane with rigid inclusions 
using Monte Carlo simulations. The phase behavior 
strongly depends on the inclusion size. For inclusion sizes 
ranging from (5 = ltoQ = 4x4 lattice sites or mem- 
brane patches, the critical stiffness K* decreases with the 
size as K* ^ Q~°'^° in the absence of an external poten- 
tial, see Fig. |g. The lateral extension of a membrane 
patch, the lattice spacing a, corresponds to the cut-off 
length for membrane shape fluctuations, which has been 
estimated as 6 nm for a membrane with a thickness of 
about 4 nm js^. In biological or biomimetic membranes, 
rigid objects with an extension larger than 6 nm may 
correspond to large trans-membrane proteins, aggregates 
of proteins and other macromolecules, or, more general. 



membrane domains with increased elastic moduli. Col- 
loidal particles adsorbed on membranes suppress mem- 
brane fluctuations similar to rigid inclusions. In general, 
membrane inclusions and membrane-adsorbed particles 
may have a variety of shapes and, therefore, orientational 
degrees of freedom |2^, |2^, |2^. Here, I have only con- 
sidered quadratic objects on a square lattice. At high 
area concentrations a^X > 0.8, the phase behavior of 
the membrane with quadratic inclusions of size Q > 1 is 
complicated by the packing transitions of the hard-square 
lattice gas |38, 3^ which are not considered here. These 
transitions are induced by the hard-square interactions 
(^ of the inclusions, but do not depend on the inclusion 
stiffness in contrast to the dynamic phase separation. 

The inclusions considered here suppress fluctuations of 
the local curvatures Ci^x and Ci^y in x- and y-direction 
at the inclusion sites, see Eqs. (3) to (5). In con- 
trast, inclusions with increased bending rigidity studied 
in ||2^, ^ only suppress fluctuations of the total cur- 
vature Ci^x + Ci^y, but not 'saddle- type' fluctuations with 
Ci,x = ~Ci,y, which seems somewhat less realistic. The 
phase behavior of these inclusions is remarkably different 
from that of the rigid inclusions studied here. Inclusions 
with increased bending rigidity do not interact in the 
absence of an external membrane potential, since local 
fluctuations of the total curvature at different membrane 
sites are not correlated in the free membrane [Q. Such 
correlations are only induced by the external potential 
of the membrane. The fluctuation-induced interactions 
between inclusions with increased bending rigidity attain 
a maximum at a certain nonzero potential strength, but 
are always considerably weaker than those of rigid inclu- 
sions characterized by Eq. (5), see Ref. ||2^. Intermediate 
cases between these two types of inclusions can be stud- 
ied by using an elastic energy with two moduli for the 
inclusions ||29|. 

Rigid inclusions may also be subject to other 
membrane-mediated interactions if they perturb the bi- 
layer thickness or have a conical or wedge-like shape, see 
introduction. In general, dynamic fluctuation-induced 
interactions can be assumed to be additive to static in- 
teractions arising from perturbations of the equilibrium 
membrane structure, as long as these perturbations do 
not affect the elastic moduli of the membrane. The 
variety of membrane-mediated indirect interactions of- 
ten complicates the interpretation of experimental re- 
sults. The dynamic phase separation of a fluctuating 
multi-component membrane in contact with a substrate 
has been recently reported in po| . The membrane con- 
tains anchored PEG-polymers and appears to phase sepa- 
rate into domains with different separation from the sub- 
strate, which might result from different effective bending 
rigidities for the domains since the fluctuation-induced 
Helfrich repulsion between membrane and substrate 
depends on the rigidity. Fluctuation-induced interac- 
tions may also contribute to the aggregation of latex 
spheres adsorbed to vesicles reported in Ref. |4^ . Lateral 
phase separation has also been observed during the ad- 
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hesion of biomimetic membranes with specific receptors 
or stickers which bind to hgands in a supported mem- 
brane 1^. Phase separation during membrane adhesion 
may be induced by membrane fluctuations |32| |33| , or 
by an effective barrier in the interaction energy between 
the membranes [Ssl 0, 44 1. In the first case, the aggre- 



gation tendency of the stickers strongly increases with 
the sticker size, similar to the rigid inclusion considered 
here The fluctuation-induced interactions between 
bound stickers are also enhanced if the stickers are more 
rigid than the surrounding membrane |33]. 
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